/* "Efficiency and water use: Dynamic effects of irrigation technology adoption"
by Micah Cameron-Harp and Nathan Hendricks

Code written by Micah Cameron-Harp
May 16th, 2024

This do file creates Figures A1 and A2 displaying the sizes of adoption cohorts.
*/

*Define directories and set working directory
/* NOTE - To replicate our results, you need to change the root directory address in the next line */
global dr_root = "C:\Users\Micah\Dropbox\Irrigation technology transition\final revisions for conditional acceptance\replication materials"
global dr_code = "${dr_root}\code"
global dr_data = "${dr_root}\data"
global dr_output = "${dr_root}\outputs"
global dr_output_main = "${dr_root}\outputs\main_text"
global dr_output_app = "${dr_root}\outputs\appendices"
global dr_output_log = "${dr_root}\outputs\logs"
global dr_temp = "${dr_root}\data\intermediate"
cd "${dr_root}"

*********************** Flood to CP or LEPA ************************************
*Load in the data for the transition from flood to cp/lepa
import delimited using "$dr_temp\flood_cporlepa_prepped.csv", clear
keep if num_switch_flood_cplepa==1
keep if wua_year <= 2015

*Create three indicators for never-treated and not-yet-treated
*'treat_status' indicates if a unit is currently treated
gen nyt = cohort_flood_cplepa > 0 & wua_year < cohort_flood_cplepa
gen nev_treat = cohort_flood_cplepa==0

*Process data
collapse (sum) treat_status nyt nev_treat, by(wua_year)

foreach var of varlist treat_status nyt nev_treat {
	gen share_`var' = 100*`var'/(treat_status+nyt+nev_treat)
}
rename (treat_status nyt nev_treat share_treat_status share_nyt share_nev_treat) ///
	(treatcat1 treatcat2 treatcat3 share_treatcat1 share_treatcat2 share_treatcat3)
	
reshape long treatcat share_treatcat, i(wua_year) j(treat_cat)
rename treatcat num_wrgs

gen treat_name = "Treated"
replace treat_name = "Not-yet-treated" if treat_cat==2
replace treat_name = "Never-treated" if treat_cat==3

label define treat_codes 1 "Treated" 2 "Not-yet-treated" 3 "Never-treated"
label values treat_cat treat_codes

*Xtset and graph
xtset treat_cat wua_year
	xtline num_wrgs, overlay
		
*Stack the data 
gen stack_wrgs  = .   // generate empty variables
gen stack_share_wrgs = .
sort wua_year treat_cat  // it is import to sort the data here.
levelsof wua_year, local(years) 
	foreach y of local years {
	summ treat_cat
	// take the first tech as it is
	replace stack_wrgs = num_wrgs if wua_year==`y' & treat_cat==`r(min)'
	// iteratively add up the 
	replace stack_wrgs = num_wrgs + stack_wrgs[_n-1]  if  wua_year==`y' & treat_cat!=`r(min)'  
	  
	* stack shares  
	replace stack_share_wrgs = share_treatcat if  wua_year==`y' & treat_cat==`r(min)' 
	replace stack_share_wrgs = share_treatcat + stack_share_wrgs[_n-1] if  wua_year==`y' & treat_cat!=`r(min)'
}

*Graph stacked share of acres
xtline stack_share_wrgs, overlay 

*Now fill in the regions and fancy graph 
**** get the labels in order
summ wua_year 
gen tag = 1 if wua_year==`r(max)'
sort wua_year treat_cat
*Generate rank variable 

*Reverse ordering 
sort wua_year treat_cat
summ treat_cat  
gen treat_cat2 = `r(max)' + 1 - treat_cat  // reverse the ordering
sort wua_year treat_cat
gen labely = .
 replace labely =  stack_share_wrgs / 2           if treat_cat==1 & tag==1
 replace labely = (stack_share_wrgs + stack_share_wrgs[_n-1]) / 2  if treat_cat!=1 & tag==1
egen rank = rank(share_treatcat)   if tag==1 & share_treatcat > 0 & share_treatcat!=., f 
 
format share_treatcat %9.1f
gen labelval = treat_name + " (" + string(share_treatcat, "%9.1fc") + "%)"   if rank <= 10
levelsof treat_cat2, local(levels)
local items = `r(r)'
 
foreach x of local levels {
display "`x'"
local x2 = `x' + 2
colorpalette mrc, n(7) nograph
 local stacklines2 `stacklines2' area stack_share_wrgs wua_year if treat_cat2 == `x', fcolor("`r(p`x2')'") lcolor(black) lwidth(*0.2) ||
 }
summ wua_year 
local x1 = `r(min)'
local x2 = `r(max)'
summ wua_year
summ treat_cat if wua_year==`r(max)'
  local ctot = `r(sum)'
  local ctot : di %7.0fc `ctot'
  
graph twoway `stacklines2' ///
 (scatter labely wua_year if tag==1, ms(smcircle) msize(0.2) mcolor(black%0) mlabel(labelval) mlabsize(vsmall) mlabcolor(black) ) ///
 , ///
  legend(off) ///
  ytitle("", size(vsmall) angle(vertical)) ///
  ylabel(, labsize(vsmall) nogrid) ///
  xtitle("") ///
  xlabel(`x1'(2)`x2', labsize(vsmall) angle(vertical) nogrid) ///
  xscale(range(1991 2021)) ///
  title("Flood to center pivot", size(small) position(11)) ///
  graphregion(margin(1 1 1 1)) xsize(4.5)
  
*Save
cd "${dr_output_app}"
graph save treat_comp_over_time_floodtocplepa, replace

*********************** CP to LEPA *********************************************
*Load in the data for the transition cp to lepa
import delimited using "$dr_temp\cp_lepa_prepped.csv", clear

*Create three indicators for never-treated and not-yet-treated
*'treat_status' indicates if a unit is currently treated
gen nyt = cohort_cp_lepa > 0 & wua_year < cohort_cp_lepa
gen nev_treat = cohort_cp_lepa==0

*Process data
collapse (sum) treat_status nyt nev_treat, by(wua_year)

foreach var of varlist treat_status nyt nev_treat {
	gen share_`var' = 100*`var'/(treat_status+nyt+nev_treat)
}
rename (treat_status nyt nev_treat share_treat_status share_nyt share_nev_treat) ///
	(treatcat1 treatcat2 treatcat3 share_treatcat1 share_treatcat2 share_treatcat3)
	
reshape long treatcat share_treatcat, i(wua_year) j(treat_cat)
rename treatcat num_wrgs

gen treat_name = "Treated"
replace treat_name = "Not-yet-treated" if treat_cat==2
replace treat_name = "Never-treated" if treat_cat==3

label define treat_codes 1 "Treated" 2 "Not-yet-treated" 3 "Never-treated"
label values treat_cat treat_codes

*Xtset and graph
xtset treat_cat wua_year
	xtline num_wrgs, overlay
		
*Stack the data 
gen stack_wrgs  = .   // generate empty variables
gen stack_share_wrgs = .
sort wua_year treat_cat  // it is import to sort the data here.
levelsof wua_year, local(years) 
	foreach y of local years {
	summ treat_cat
	// take the first tech as it is
	replace stack_wrgs = num_wrgs if wua_year==`y' & treat_cat==`r(min)'
	// iteratively add up the 
	replace stack_wrgs = num_wrgs + stack_wrgs[_n-1]  if  wua_year==`y' & treat_cat!=`r(min)'  
	  
	* stack shares  
	replace stack_share_wrgs = share_treatcat if  wua_year==`y' & treat_cat==`r(min)' 
	replace stack_share_wrgs = share_treatcat + stack_share_wrgs[_n-1] if  wua_year==`y' & treat_cat!=`r(min)'
}

*Graph stacked share of acres
xtline stack_share_wrgs, overlay 

*Now fill in the regions and fancy graph 
**** get the labels in order
summ wua_year 
gen tag = 1 if wua_year==`r(max)'
sort wua_year treat_cat
*Generate rank variable 

	*Reverse ordering 
	sort wua_year treat_cat
	summ treat_cat  
	gen treat_cat2 = `r(max)' + 1 - treat_cat  // reverse the ordering
sort wua_year treat_cat
gen labely = .
 replace labely =  stack_share_wrgs / 2           if treat_cat==1 & tag==1
 replace labely = (stack_share_wrgs + stack_share_wrgs[_n-1]) / 2  if treat_cat!=1 & tag==1
egen rank = rank(share_treatcat)   if tag==1 & share_treatcat > 0 & share_treatcat!=., f 
replace rank = 2 if tag==1 & treat_cat2==2
replace rank = 3 if tag==1 & treat_cat2==1
*Make space for each label
replace labely = 94.5 if rank==2
replace labely = 99 if rank==3


 
format share_treatcat %9.1f
gen labelval = treat_name + " (" + string(share_treatcat, "%9.1fc") + "%)"   if rank <= 10
levelsof treat_cat2, local(levels)
local items = `r(r)'
 
foreach x of local levels {
display "`x'"
local x2 = `x' + 2
colorpalette mrc, n(7) nograph
 local stacklines2 `stacklines2' area stack_share_wrgs wua_year if treat_cat2 == `x', fcolor("`r(p`x2')'") lcolor(black) lwidth(*0.2) ||
 }
summ wua_year 
local x1 = `r(min)'
local x2 = `r(max)'
summ wua_year
summ treat_cat if wua_year==`r(max)'
  local ctot = `r(sum)'
  local ctot : di %7.0fc `ctot'
  
graph twoway `stacklines2' ///
 (scatter labely wua_year if tag==1, ms(smcircle) msize(0.2) mcolor(black%0) mlabel(labelval) mlabsize(vsmall) mlabcolor(black) ) ///
 , ///
  legend(off) ///
  ytitle("", size(vsmall) angle(vertical)) ///
  ylabel(, labsize(vsmall) nogrid) ///
  xtitle("") ///
  xlabel(`x1'(2)`x2', labsize(vsmall) angle(vertical) nogrid) ///
  xscale(range(1991 2025)) ///
  title("Traditional center pivot to LEPA", size(small) position(11)) ///
  graphregion(margin(1 1 1 1)) xsize(4.5)
 
*Save
cd "${dr_output_app}"
graph save treat_comp_over_time_cptolepa, replace

**************** Combine graphs into Figure A3 *********************************
*Bring in the treatment composition graphs 
graph use "${dr_output_app}\treat_comp_over_time_floodtocplepa.gph", name(flood_treat_comp, replace)
graph use "${dr_output_app}\treat_comp_over_time_cptolepa.gph", name(cplepa_treat_comp, replace)

*Combine them
graph combine flood_treat_comp cplepa_treat_comp,  ///
	cols(1) name(comb_gr, replace) ///
	title("", size(medsmall)) ///
	b2title("Year", size(vsmall)) l2title("Percent of sample", size(vsmall)) ///
	xsize(4.5) ysize(4.5)
	
graph export "${dr_output_app}\figureA3.tif", replace wid(6500) height(6500)